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Abstract. 

Theoretical models suggest that the first stars in the universe could have been very 
massive, with typical masses > 100 M©. Many of them might have died as energetic 
thermonuclear explosions known as pair-instability supernovae (PSNe). We present 
multidimensional numerical simulations of PSNe with the new radiation-hydrodynamics 
code CASTRO. Our models capture all explosive burning and follow the explosion un- 
til the shock breaks out from the stellar surface. We find that fluid instabilities driven by 
oxygen and helium burning arise at the upper and lower boundaries of the oxygen shell 
~ 20 - 100 sec after the explosion begins. Later, when the shock reaches the hydrogen 
envelope a strong reverse shock forms that rapidly develops additional Rayleigh-Taylor 
instabilities. In red supergiant progenitors, the amplitudes of these instabilities are suf- 
ficient to mix the supernova's ejecta and alter its observational signature. Our results 
provide useful predictions for the detection of PSNe by forthcoming telescopes. 



1. Introduction 



The evolution of the first stars in the universe is one of the frontiers of modem cos- 
mology. Primordial stars synthesized the first heavy elements in the universe, and their 
energetic feedback influenced th e formation of la ter generations of stars and the first 
galaxies (^ Whalen et all 12008 a. 'b^: iGreif et alj |201(/). Ea rly numerical models predicted 
that P op III stars formed with masses of 100-1000 Mq dBromm et al.ll2009. : Abel et all 
2002h. Ne\y studies have found that ~ 20% of Pop III stars form in binaries or multiples 
( Turk et aLll2009l : IStacy et aLll2010l) so the first stars could be less massive than orig- 
inally thought. However, e ven today observatio ns support the existence of stars with 
initial masses over 150 M© ( Crowther et al.ll2010h . Stellar evolution models predict that 
Pop III stars with initial masses of 140 - 260 Mp d evelop oxygen cores of > 50 Mq after 
central carbon burning (iHeger & Woosleyll200 2l). At this point the core reaches suffi- 
ciently high temperatures (~ 10^ K) and low densities (~ 10^ g/cc) that the creation of 
electron-positron pairs is favored. Radiation pressure support then quickly decreases, 
triggering a rapid contraction of the core. During contraction, core temperatures and 
densities sharply rise and oxygen and silicon begin to burn explosively. The result- 
ing thermonuclear explosion, known as a pair-instability supernova (PSN), reverses the 
contraction and completely unbinds the star, leaving no compact remnant and forming 



1 



2 



Chen, Heger, and Almgren 



up to 50 of^^Ni. One possible PSN candidate, SN 2007bi, has recendy found by 
(lGal-YametalJl2009h . 

Most current theoretical models of PSNe are based on one-dimensional calcula- 
tions ( Heger & Woosleyl 12002 ). However, in the initial stages of a supernova spherical 
symmetry is broken by fluid instabilities generated by burning, which cannot be cap- 
tur ed in ID. Two-dim e nsion al simulations of Pop III PSNe have recently been done 
by lJoggerst & Whaled (l201ll) in which only mild dynamical instabilities were found to 
form, but they proceeded from ID KEPLER models in which explosive burning had 
already occurred and thus exclude instabilities driven by burning. Such instabilities, if 
they form, may alter the energetics and nucleosynthesis of the SN by vigorously mixing 
its fuel and must be included in simulations to understand the true evolution of PSNe. 
We have performed 2D simulations of Pop III PSNe that follow the initial contraction 
of the core until most of the energy due to explosive burning has been released, in 
contrast to IJoggerst & WhalenI (l20Ilb . who only follow the post-nucleosynthesis hy- 
drodynamics. Our goal is to study any fluid instabilities that arise and how mixing 
alters nucleosynthesis and the energetics of the explosion. 



2. Numerical Methods 



We evolve zero-metallicity stars in KEPLER (iWeaver et al.lll978h . a one-dimensional 



Lagrangian stellar evolution code. In KEPLER we solve evolution equations for mass, 
momentum, and energy and include physics relevant to stellar evolution such as nuclear 
burning and artificial mixing. When the star comes to the end of central oxygen burning, 
we map its profile onto a 2D Cartesian grid in CASTRO. The procedure f or mapping 



and see ding initial perturbations in these profiles is discussed in detail in iChen et al. 
(1201 la ). We evolve the star in CASTRO through the end of explosive burning. 



CASTRO (Almgren et al. 2010 ; Zhang et al. 2011) is a massively parallel, multidi- 
mensional Eulerian adaptive mesh refinement (AMR) radiation-hydrodynamics code 
for astrophysical applications. Its time integration of the hydrodynamics equations is 
based on a higher-order, unsplit Godunov scheme. Block-structured AMR with subcy- 
cling in time enables the use of high spat ial resolution where it is most needed. We use 
the Helmholtz equation of state (EOS) (ITimmes & Swe stv"200(f) with density, temper- 
ature, and species mass fractions as inputs. The gravitational field is calculated using a 
monopole approximation constructed from a radial average of the 2D density field on 
the grid. 



3. 2D Simulations 

In Fig. [T] we show the formation of dynamical instabihties at the base of the oxygen 
burning shell during the contraction of the core dChen et al.ll2011bl) . They are relatively 
mild and do not penetrate the central ^^Ni region, so no ^^Ni is mixed into the upper 
layers of the star at this stage. After explosive burning reverses the contraction of the 
core, fluid instabilities driven by helium burning also appear in the outer layers of the 
oxygen shell. Minor mixing caused by these instabilities begins about 100 sec after 
reversal of the collapse. 

As we show in Fig. |2j when the shock propagates into the hydrogen envelope 
the formation of a strong reverse shock creates additional Rayleigh-Taylor instabilities 
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Figure 1. Fluid instabilities at early stages of the explosion. We show oxygen 
mass fraction 100 sec after the reversal of core collapse. Dynamical instabilities ap- 
pear at the lower boundary of the O shell because of oxygen burning and at the upper 
boundary of the shell because of He burning. The red lines indicate the location of 
shock. 

(RTI). Their amplitudes are sufficient to mix oxygen with the surrounding shells: H, 
He, and Si. Some mixing also occurs at the outer edge of the ^^Ni core. Our results 
demonstrate that dynamical instabilities form at several stages of the explosion and that 
they are mainly driven by RTI at the interfaces of contact discontinuities in density or 
species abundance. Mixing is important to the observational signatures of Pop III PSNe 
because it can cause absorption lines of heavy elements to appear in spectra sooner 
after shock breakout than when mixing is absent. The instabilities can also lead to the 
formation of clumps in the ejecta that can strongly affect its luminosity at later times. 

4. Conclusions 

In contrast to earlier multidimensional simulations, we find that mixing can occur at 
several stages in Pop III PSNe prior to shock breakout when core contraction and ex- 
plosive nuclear burning is followed in 2D. Mixing appears in PSN spectra by introduc- 
ing absorption lines of heavy elements at early times and by changing its luminosity. 
Our 2D simulations are the first in a numerical campaign to investigate the evolution of 
PSNe from their earliest stages as pre-supernova progenitors. We are currently prepar- 
ing to post-process our simulations to compute light curves and spectra in order to pro- 
vide useful predictions for the James Webb Space Telescope, which may soon detect 
these primordial explosions. 
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Figure 2. Fluid instabilities prior to shock breakout. Here, we show the density 
when the shock is about to break out from the stellar surface. The fluid instabilities 
are driven by the reverse shock and lead to significant mixing. 
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